home *** CD-ROM | disk | FTP | other *** search
/ Practical Algorithms for Image Analysis / Practical Algorithms for Image Analysis.iso / TARFILE.GZ / tarfile / ch_6.2 / xptstats / xptstats.c < prev    next >
Encoding:
C/C++ Source or Header  |  1999-09-11  |  6.0 KB  |  259 lines
  1. /* 
  2.  * xptstats.c
  3.  * 
  4.  * Practical Algorithms for Image Analysis
  5.  * 
  6.  * Copyright (c) 1997, 1998, 1999 MLMSoftwareGroup, LLC
  7.  */
  8.  
  9. /*
  10.  * XPTSTATS
  11.  *
  12.  * Asses the degree of randomness in a point pattern
  13.  * by drawing and filling circles of increasing radius
  14.  * whose centers lie at the specified points.
  15.  *
  16.  * input: 
  17.  * A .vin file (generated by spp)
  18.  *
  19.  * output:
  20.  * a
  21.  */
  22.  
  23. #include "xptstats.h"
  24.  
  25. /* globals */
  26. extern char *optarg;
  27. extern int optind, opterr;
  28. extern short tiffInput;         /* flag=0 if no ImageIn to set tags; else =1 */
  29.  
  30. /*
  31.  * readpoints()
  32.  *   DESCRIPTION:
  33.  *     read points from file
  34.  *   ARGUMENTS:
  35.  *     filename: string for file to read from
  36.  *   RETURN VALUE:
  37.  *     none
  38.  */
  39.  
  40. void
  41. readpoints (char *filename, int *npoints, struct Point **points, float *xmin, float *xmax, float *ymin, float *ymax)
  42. {
  43.   int i;
  44.   char buf[1024];
  45.   char *strptr;
  46.   FILE *file;
  47.  
  48.   strcpy (buf, filename);
  49.   if (((strptr = strstr (buf, FILE_EXT)) && (strptr != (buf + strlen (buf) - strlen (FILE_EXT)))) || !strptr)
  50.     strcat (buf, FILE_EXT);
  51.   if ((file = fopen (buf, "r")) == NULL) {
  52.     printf ("couldn't open file %s!!\n", buf);
  53.     exit (1);
  54.   }
  55. /*
  56.  * skip top line in .vin data file
  57.  */
  58.   fscanf (file, "%d %f %f %f %f", npoints, xmin, xmax, ymin, ymax);
  59.  
  60. /*
  61.  * read data
  62.  */
  63.   *points = (struct Point *) calloc (*npoints, sizeof (struct Point));
  64.   for (i = 0; i < *npoints; i++)
  65.     fscanf (file, "%f %f", &(*points + i)->x, &(*points + i)->y);
  66.  
  67.  
  68. #ifdef SHOW_DATA
  69.   printf ("\n...points read in are:\n");
  70.   printf ("...parameters: npoints: %4d", *npoints);
  71.   printf ("\n   extrema: %f %f %f %f\n", *xmin, *xmax, *ymin, *ymax);
  72.  
  73.   for (i = 0; i < *npoints; i += 1) {
  74.     printf ("  x[%d] = %f   y[%d] = %f\n",
  75.             i, (*points + i)->x, i, (*points + i)->y);
  76.   }
  77. #endif
  78.   fclose (file);
  79. }
  80.  
  81.  
  82. void
  83. exit_messg (str, code)
  84.      char *str;
  85.      int code;
  86. {
  87.   printf ("\n%s", str);
  88.   exit (code);
  89. }
  90.  
  91. /*
  92.  * usage of routine
  93.  */
  94. void
  95. usage (char *progname)
  96. {
  97.   progname = last_bs (progname);
  98.   printf ("USAGE: %s infile [-t t_max] [-n ns_max] [-w file] [-L]\n", progname);
  99.   printf ("\n%s constructs the distribution functions N(t)/N_0\n");
  100.   printf ("and P(t) for a given point pattern\n\n");
  101.   printf ("ARGUMENTS:\n");
  102.   printf ("        infile: input filename (.vin) (ASCII)\n");
  103.   printf ("                produced by spp or xah\n\n");
  104.   printf ("OPTIONS:\n");
  105.   printf ("   -t t_max: maximum radius (R/R0) at which to stop (default = %d)\n", T_MAX_DEF);
  106.   printf ("  -n ns_max: maximum # iteration steps (default = %d)\n", NS_MAX_DEF);
  107.   printf ("    -w file: write data to file\n");
  108.   printf ("         -L: print Software License for this module\n");
  109.   exit (1);
  110. }
  111.  
  112.  
  113. void
  114. main (int argc, char *argv[])
  115. {
  116.   int i;
  117.   int i_arg;
  118.   int ns;
  119.   int ns_max = NS_MAX_DEF;
  120.   int nx;
  121.   int ny;
  122.   int S0;
  123.   int N0;
  124.   float t;
  125.   float t_max = (float) T_MAX_DEF;
  126.   float xmin;
  127.   float xmax;
  128.   float ymin;
  129.   float ymax;
  130.   float R;
  131.   float R0;
  132.   float DelR;
  133.   float *t_arr;
  134.   float *ncr_arr;
  135.   float *p_arr;
  136.   int fill_index;
  137.   int circle_index;
  138.   int n;
  139.   long n_covered;
  140.   struct Point *points;
  141.   int npoints;
  142.   Image *imgOut;
  143.   char outFile[256];
  144.   int write_flag = OFF;
  145.   FILE *fpOut;
  146.  
  147. /* 
  148.  * cmd line options
  149.  */
  150.   static char *optstring = "t:n:w:L";
  151.  
  152. /*
  153.  * parse command line
  154.  */
  155.   optind = 2;
  156.   opterr = ON;                  /* give error messages */
  157.  
  158.   if (argc < 2)
  159.     usage (argv[0]);
  160.  
  161.   while ((i_arg = getopt (argc, argv, optstring)) != EOF) {
  162.     switch (i_arg) {
  163.     case 't':
  164.       sscanf (optarg, "%f", &t_max);
  165.       break;
  166.     case 'n':
  167.       sscanf (optarg, "%d", &ns_max);
  168.       break;
  169.     case 'L':
  170.       print_sos_lic ();
  171.       exit (0);
  172.     case 'w':
  173.       strcpy (outFile, optarg);
  174.       write_flag = ON;
  175.       break;
  176.     default:
  177.       printf ("...unknown command line arg encountered\n");
  178.       exit (1);
  179.       break;
  180.     }
  181.   }
  182.   /* reset tiffInput so that we write a grayscale file (i.e tags are not copied) */
  183.   tiffInput = 0;
  184.  
  185.   readpoints (argv[1], &npoints, &points, &xmin, &xmax, &ymin, &ymax);
  186.  
  187.   S0 = (int) ((xmax - xmin) * (ymax - ymin));
  188.   N0 = npoints;
  189.   fill_index = 127;
  190.   circle_index = 255;
  191.   R0 = (float) (sqrt ((float) S0 / (float) N0));
  192.   R = (float) 0.0;
  193.   DelR = R0 / (float) ns_max;
  194.  
  195.   printf ("Parameters are:\n");
  196.   printf ("t_max = %f\n", t_max);
  197.   printf ("S0 = %d, N0 = %d, R0 = %f, ns_max = %d, DelR = %f\n", S0, N0, R0, ns_max, DelR);
  198.  
  199.  
  200.   t_arr = (float *) calloc (ns_max, sizeof (float));
  201.   ncr_arr = (float *) calloc (ns_max, sizeof (float));
  202.   p_arr = (float *) calloc (ns_max, sizeof (float));
  203.  
  204. /*
  205.  * initialize graphics
  206.  * (allocate enough image area for drawing filled circles!!)
  207.  */
  208.   imgOut = ImageAlloc ((long) (ymax + ymin), (long) (xmax + xmin), BPS);
  209.  
  210.   ns = 0;
  211.   while ((t = R / R0) < t_max) {
  212.     reset_image (imgOut, 0);
  213.     R += DelR;
  214.     for (i = 0; i < N0; i++) {
  215.       draw_filled_circle ((int) points[i].x, (int) points[i].y, (int) R, imgOut, circle_index);
  216.     }
  217.  
  218.     n = 0;
  219.     n_covered = 0;
  220.     for (i = 0; i < N0; i++) {
  221.       if (getpixel ((int) points[i].x, (int) points[i].y, imgOut) == circle_index) {
  222.         gdImageFill ((int) points[i].x, (int) points[i].y, imgOut, fill_index);
  223.         n++;
  224.       }
  225.     }
  226.     n_covered = 0;
  227.     for (nx = (int) xmin; nx < (int) xmax; nx++) {
  228.       for (ny = (int) ymin; ny < (int) ymax; ny++) {
  229.         if (getpixel (nx, ny, imgOut))
  230.           n_covered++;
  231.       }
  232.     }
  233.     t_arr[ns] = R / R0;
  234.     ncr_arr[ns] = (float) n / (float) N0;
  235.     p_arr[ns] = (float) n_covered / (float) S0;
  236.     if (ns++ >= ns_max)
  237.       break;
  238.   }
  239.  
  240.  
  241.   printf ("\nThe distribution functions are:\n\n");
  242.   printf ("   T      N(T)/N0    P(T)\n");
  243.   for (i = 0; i < ns; i++)
  244.     printf ("%f %f %f\n", t_arr[i], ncr_arr[i], p_arr[i]);
  245.  
  246.   if (write_flag) {
  247.     if ((fpOut = fopen (outFile, "w")) == NULL) {
  248.       printf ("\n...cannot open %s for writing\n", outFile);
  249.       exit (1);
  250.     }
  251.     printf ("\n...writing t, N(t)/N_0 and P(t) to %s\n", outFile);
  252.     for (i = 0; i < ns; i++)
  253.       fprintf (fpOut, "%f %f %f\n", t_arr[i], ncr_arr[i], p_arr[i]);
  254.     fclose (fpOut);
  255.   }
  256.   ImageFree (imgOut);
  257.  
  258. }
  259.